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ABSTRACT 

A Liquid Metal Fast Breeder Reactor is simulated on a hybrid 
canputer to study the transient behavior of the neutron density during 
a controlled reactivity input disturbance. 

The nonlinear partial differential equations of heat flow are 
reduced by a discrete space-continuous time method to ordinary non- 
linear differential equations v\^ich are readily solved on the analog 
canputer. Use is made of time multiplexing of the analog circuitry 
in order to reduce the number of conponents. An open-loop iteration 
process is employed to solve the closed- loop feedback system. 

Recently conducted research with the open- loop iteration method 
has demonstrated that a large number of iterations are required for 
convergence. An algorithm is developed which gives an improvement in 
the convergence rate for early values of time. For times greater than 
two seconds a stability problem with the converged solution was en- 
countered and is discussed with seme observations and comments. 

Innovations are introduced in the simulation of the neutron kinetics 
equations and in the handling of the nonlinear thermal properties of the 



uranium oxide fuel. 
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Reactivity units 
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Fuel expansion coefficient 


Reactivity units/^C 
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Sodium void coefficient 


Reactivity units/°C 


c^(t) 


th 

Density of the i group of 
delayed neutron precursor 


3 

Number/cm 


c(t) 


Density of the assumed one 
precursor at time t 


Number/cm^ 


c(0) 


Density of the assumed one 
precursor at time t = 0 


3 

Number/cm 


c 


Specific heat (subscript f 
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and cool to coolant) 
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d 


Density (subscript f refers 
to fuel, c to cladding, and 
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H 
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Symbol Definition Units 

Ar Radial width of node within Cm 

fuel 



Ar Radial width of cladding Qn 

c 

Ar , Radial width of coolant Cm 

around fuel pin 



R 



s(z) 

S(R) 



The three spacial variables Cm 

of core 

3 
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3 

Spacial power distribution Watts/on 
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t 
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1 
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T 
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c 



T 
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-1 
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AT 
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AT^ 



C 



Time Seconds 

Temperature (subscript i 
refers to radial node and 
superscript j refers to 
axial node) 

Temperature of radial node i *^K 
and axial node j at operating 
level 

Av^age fuel temperature of 
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level 
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Change in average coolant °K 

temperature of axial 

segment 



11 



Symbol 


Definition 


Units 
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tanperature of axial 
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AT 
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°K 


4 
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4 


Coolant weighting factor of 
jth axial segment 
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Z 


Axial variable of fuel pin 


Cm 


Az 


Height of axial node 


Cm 




Delayed neutron fraction from 
i^ precursor 


Pure number 


B 


Delayed neutron fraction 


Pure number 


X. 

1 


Decay constant for i^ 
precursor 


Sec ^ 


A 


Decay constant for the assumed 
one precursor 


Sec ^ 


Pe 


Reactivity due to controlled 
disturbance 


Reactivity units 


"f 


Reactivity due to feedback 


Reactivity units 


P 


Reactivity (subscript Dopp 
refers to Doppler effect, 
exp to fuel expansion, and 
cool refers to coolant 
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Reactivity units 


ip(R) 
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density 
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I . INTRDDUCTION 



For over twenty years the nuclear community in the U. S. has been 
working to meet the recognized important objectives of econcmic nuclear 
power and improved fuel utilization through the developnent of the 
Liquid Metal Fast Breeder Reactor [Ref. 1]. This effort, however, 
has been largely dispersed with no real organized effort. Starting 
in 1965 an unprecedented demand for light water reactor power plants 
occurred paralleled with increased loranium needs. This increased fuel 
requirement emphasized the need to develop a Liquid Metal Fast Breeder 
Reactor and has caused the U. S. as well as other countries to establish 
Liquid Metal Fast Breeder Reactor (LMFBR) programs. 

The light water reactors are very successful today with favorable 
prediction of continued improvements in fuel performance and reduced 
operating cost [Ref. 2] . However, the IMFBR can provide the most ef- 
ficient means of exploiting the energy available in uranium because of 
a high thermodynamic efficiency, an efficient fuel conversion cycle, 
reasonable capital and operating costs, and safety in operation. 

The present U. S. effort in Fast Reactors is summed up in the 
nuclear reactors Clementine (0.025 Megawatts), LAMPRE I (1 Megawatt), 
EBR-I (1.4 Megawatts), EBR-II (62.5 Megawatts) and Enrico Fermi (200 
Megawatts) of vdiich only the latter two remain in canmission. The 
current LMFBR program envisions a 1,000 Megawatt electric reactor, 
v^ich is a substantial increase over the largest systen constructed 
to date. It is generally accepted [Ref. 3 and 4] that the operating 
condition of the proposed LMFBR is beyond the range vdiich permit prudent 
extrapolation from existing experiences. Simulation studies can greatly 
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assist in design evaluation by predicting the behavior of the conposite 
model based on experimental information obtained from the ccnponent 
parts. 

The use of hybrid simulation is a very effective method to carry 
out design evaluation based upon the transient response to predictable 
accidental disturbances. The pure analog approach is not practical due 
to the prohibited amount of hardware required to simulate the coupled 
nonlinear differential equations of heat flow. The pure digital ap- 
proach, which is feasible, does not readily lend itself to multiple 
parameter studies for design optimization; the computer time required 
would be economically unfeasible. The hybrid computer combines seme 
of the best characteristics of the analog and digital machines result- 
ing in a more econcmical computer. Recently, large-scale simulations 
[Ref. 5 and 6] were performed with the object of comparing the economics 
of hybrid versus digital for large plant simulation. These simulations 
indicate that this advantage is not about to disappear. 

Previous efforts have been made to daronstrate the feasibility of 
employing a hybrid computer for the analysis of the transient behavior 
of a IMFBR core [Ref. 7] . The developnent of new techniques to enhance 
the effectiveness of the hybrid computer are onarrently in progress 
[Ref. 8] . The method developed in this report is a contribution in 
this effort with regards to software design related to reducing the 
number of iterations for convergence. The number of iterations is 
critical in a hybrid simulation since a large number of iterations 
would defeat the purpose of using the hybrid computer. A pure digital 
approach would then become more practical where both discrete time and 
discrete space could be employed. 
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Additional contributions were realized by simulating the nonlinear 
thermal conductivity of UO 2 on a Diode Function Generator resulting in 
a reduction of analog hardware, and simplification of the neutron 
kinetics equations for use on the analog computer. 

The parameters for the model used are primarily from a recently 
proposed 1000 MWe LMFBR design [Ref. 9] . However the model developed 
is quite general and can be employed in the analysis and design studies 
of fast, intermediate or thermal reactors. 
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II. DESCRIPTION OF MODEL 



A. REACTOR KINETICS 

The point reactor kinetics equations are used to describe the 
neutron density behavior. The assumption is made that the density 
(F) can be separated into a time dependent and a space dependent 
term, and can be expressed as^ 



Reference 10 and 11 present good treatments of the separation of 
time and space of the neutron density. n(t) is given by the fol- 
lowing equations 



A treatment of the point reactor kinetics can also be found in Ref. 

10 and 11. 

In most applications of the above equations it is considered 
appropriate to treat N = 6, that is 6 delayed neutron precursors. 
However, in this model due to lack of sufficient analog hardware only 
one delayed precursor will be considered; its parameters will be 
appropriate averages of those of the six. 



refers to the three special variables of the core and not to 
be confused with r, the radial variable of one fuel pin. 




( 1 ) 




(2) 



where 
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Equation (2) is recast into a more convenient form by defining 
two dimensionless quantities (see Appendix A for details) 









With these definitions and the one precursor assuirption equation (2) 
becomes 






(3a) 



j; D(i) 
dt 



- X Nj (!) - b (t) ] 



(3b) 



with initial conditions N(0) = D(0) =1, and where the term K(t) 
equals and is defined in dollar units, ^(t) is the sum of 
the forcing fionction (t) , plus the feedback reactivity ^ (t) , 
or simply 



The quantity for a fast reactor is of the order 10^, thus equation 

(3a) requires a high gain integrator for similation on the analog ccrputer. 
This introduces some problems in the simulation which will be discussed 
in Section V. If equation (3a) is rearranged as follows, 

0 ^ )M7; + b(i) - {^/^)ULi) 

the term [4^) represents a small quantity compared to the other 

three terms of the equation. Neglecting this term is a reasonable ap- 
proximation and equations (3a) and (3b) become 

N(t)^ 



(a) 

(5) 

(b) 
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Equation (5b) is an initial value ordinary differenticil equation and 
its solution is easily handled by the analog ccrtputer. 

The forcing function is an cirbitrary selected but known function 
of time. In this report a terminated ramp is used to model the with- 
drawal of a control rod. The feedback reactivity is not a krOftai 
function of time but depends on the average tanperature of the reactor. 
In a fast reactor there are three main contributions to the feedback 
[Ref. 12] , the Doppler effect, fuel expansion and coolant expansion. 
The expression for feedback due to these three effects can be repre- 
sented as^, 

or, in terms of dollar units 



1 






( 6 ) 



> ^ 

, B and C , are positive constants. The temperatures 
Ttopp exp cool ^ ^ 

indicated in equation (6) represent average core temperatures. 

In this ircdel a single average fuel cell is assumed to give all 

the terrperature information necessary to construct the average core 
3 

temperature . 

The power per unit volume supplied by the fission process can be 
expressed as 

Q = Hit) SLR) 

where S (R) is proportional to If (R) . The proportionality factors include 



2 

See Appendix A for more details. 

3 

References 8, 13 and 14 all employ a single "average" fuel cell 
to compute the reactivity feedback. Reference 7 uses five fuel cells 
properly weighted to ocmpute the reactivity feedback. 
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H (0) , (average fission cross section) , E (energy released per 



fission) , and V (the average velocity of the neutrons) . In this 
model and V are treated as constants. 

B. HEAT ELCW 

Large power densities exist within the core of a 1000 MWe IMEBR 
giving rise to steep temperature gradients and extreme heat flew. 
Furthermore the thermal conductivity of the ceramic fuels being 
considered for the 1000 MWe LMFBR (namely, oxides, nitrites and 
carbides) vary nonlinear ly at the proposed high operating temperatures. 
These factors dictate that the problem of heat flow be approached on 
a microscopic basis. 

A LMFBR is composed of thousands of pencil size fuel pins approxi- 
mately three feet long in the shape of right cylinders; the coolant 
flews coaxially. A sample fuel pin illustrated in Fig. 1 is cemposed 
of a cylindrical shaped fuel assembly with a concentric shell called 
cladding. 

The heat flew equation for a cylindrical shaped body with varying 
thermal properties is; 



In a fast reactor the neutron flux is nearly constant within the fuel 
pin in the radial direction (in contrast to a thermal reactor vhere flux 
depression ocezurs within the fuel) . This faert gives rise to cylindrical 
synmetry and if axial flow of heat is neglected, v^ich is consistent 
with other models proposed [Ref. 5, 13 and 14] , equation (8) for any 
fuel pin reduces to: 




(8) 
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Figure 1. Diagram of Fuel Pin 
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(9) 



( 10 ) 



where the assumptions outlined in Section IIA for the pcv;er generation 
Q, have been introduced (i.e, the space and time separation of neutron 
density) . The term S(R) has been replaced by s(z) to indicate the 
heat generation rate in a single fuel pin is a function of the axial 
direction only. 

Since the analog ccmputer can handle only ordinary differential 
equations, equations (9) and (10) are reduced to initial value ordinary 
differential equations by employing a method of integrating over the 
spacial variable [Ref. 15] . The technique involves discretizing the 
space variable in the radial direction, ^ , into five nodes within the 
fuel and one within the c] adding (see Fig. 2) . In carrying out the 
integration the temperature is assumed to vary linearly between nodes. 

For more details of the derivations refer to Appendix A. 

As a result of employing this technique equations (9) and (10) become 



(i indicates the node which varies in number from 1-5 within the fuel) 





( 11 ) 



^ ^ ^ I = thermal conductivity at boundary between node i and 
i + 1 
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Fig\ire 2. Temperature Profile Assumed Within Fuel Pin. 
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( 12 ) 




(Superscript j refers to the axial node) 



Due to synmetry about r = 0, the condition that is applied 

as a boundary condition. Also at the clad-fuel interface the continuity 
of heat flew requires that 



K 











At the boundary between cladding and coolant the conservation of 
energy is applied to arrive at the following equation (see Appendix A 
for derivation) . 



^ (M IJc ^ 2 /V 

4o„ 

G(t)= mass flow rate of coolant 

H = heat flux across the clad-coolant interface 
r^ E radius of fuel pin (fuel + cladding) 

Again an ordinary differential equation is needed. Here, however, 
the use of finite differencing of the space variable is more practical. 
This simply requires replacing by V "T j/k? and assuming 

the coolant temperature varies linearly over length AZ. Thus equation 
( 13) can be expressed as : 

“il- /j-i H 

~I2 = average coolant temperature of j axial segment = - ^ 
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At the clad-coolant boundary the heat flow per unit time (H) 
is represented by h ~ ) vhich also equals ) - 

h is the heat transfer coefficient between clad and coolant. 

After incorporating the boundary conditions and the physical 
constants listed in Appendix C the following coupled ordinary dif- 
ferential equations are obtained. 

i [T,‘., 2 o 54 T’ ..IZSi-rJ. 

id) f^u) 

- (y-V)]f, 433 NCt) 

- + .433 Ut)N(.i) 



ac) 

on 

on 



Jf LV-l-.iayoT^\.4K-Z l‘]^ 

- [Ki.it CV-V)] i- .sia i(z)^a) 

Jj [T; + l.35;,-T,' f.SSi-Tj].- 1741 [V-i‘] 



iy - - lft,3Ti' + IIB.sX + U.sX'‘ 
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The inlet temperature to the first axial segment is knovn and 
treated as a constant during any transient. In nuclear reactors 
the inlet temperatures are usually kept constant, and any rapid 
pcwer change will have leveled off before any change occurs in the 
inlet temperature. For each subsequent segments the inlet tanper- 
ature is a known function of time obtained from analysis of the 
previous section (see Digital Software section for details) . The 
steady state special pcwer distribution s(z) is known a priori, 
therefore, equations (15) through (21) can be solved on the analog 
ccraputer . 

The bonding between the fuel and cladding is ignored in this 
model and the fuel is assumed to be in constant contact with the 
cladding. In most liVIFBR under design the bonding agent is helium. 
Neglecting this bonding will produce lower than normal temperatures 
within the fuel. 
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III. HYBRID COMPUTER TECHNIQUE 



A. DESCRIPTION OF ANALOG COyiPUrER HARDWARE 

Fran the heat flow equations given in the previous section it can 
be seen that the model requires the simulation of a large number of 
ordinary differential equations. However, an important aspect to keep 
in mind is that the description of each axial segment incorporates 
the same differential equations, the coolant inlet temperature and 
the fuel heat generation rate being the only quantities varying fron 
segment to segment. With this in mind a technique of multiplexing 
that is, time shearing, of the analog circuitry can be employed. The 
details of the multiplexing used will be explained in Section IIIB. 

To use the technique, only one set of equations (15) through (21) 
need be set up on the analog conputer. In order to solve equations 
(15) through (21) on the analog computer they must first be magnitude 
scciled to insure that the outputs do not exceed the maximum analog 
conputer voltage (+ 100^) . After proper scaling equations (15) 
through (21) can be expressed as equations (41) through (47) which 
are shewn in Appendix B. 

The approach used in setting up the solution to equations (41) 
through (47) on the analog conputer can best be explained by taking 
one equation, 

4 + . u-»rT,' ]= ^3,/ fa j (T^-TJ)] 

^ - — ' _ - ^ (42) 

.2/65" 5(i) N(t) 



28 



The right side of equations (42) represents inputs to an 
integrator whose analog symbol is shown by elorvent A in Fig. 3. 
Element B is a summer amplifier. This process is extended over 
the remaining heat equations. Figure 4a is a diagram of the analog 
circuitry. 




I.C. refers to initial conditions. 

Figiare 3. Illustration of Technique Used to Construct Analog 
Circuit. 

The nonlinearity of the thermal conductivity, K(T) , was readily 
handled by the analog ccmputer through utilization of a Diode Function 
Generator (DFG) . A DFG can be used to approximate a nonlinear curve 
by the method of straight line segments. The DFG used with this 
simulation has ten break points. Thus a maximum of ten straight 
lines were available to duplicate the curve and this was quite adequate 
to give a good representation of the nonlinearity. Experimental data 
for the K(T) curve was obtained from Ref. 16 and 17. Figure 5 shows 
both curves. The DFG operates sirtply as shewn in Fig. 6. The input 
is a voltage vhich represents °K temperature and the output is a volt- 
age proportional to the thermal conductivity for that tonperature 



29 




Figure 4a. Diagram of Analog Circuit. 
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Figxire 4b. Diagram of Multiplexed DPG. 
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Figxire 4c. Diagram of Logic Circuit. 



32 




Integrator 




Summer amplifier 




Track-hold amplifier 




Analog multiplier 



o~ 



Potentiometer 




Delay flop 



Nand-gate 




Diode Function Generator 



Analog/ digital switch 



Input or output from digital computer 



Figxire 4d. Key For Analog Syitibols. 



33 




o 



o 




i3oo /7oo Zioa 2foo 

"Tayw pe val o’i'e — ®K 







VoLTs 

9oii /3od /^oo 2/00 Sj-ao 



Tevn pe'T aiure. — 



Figure 5. Thermal Conductivity Of UO^. 



34 



T _ 

(\fcii:) 



DF6 



^ K(T) 
(Volts) 



Figure 6. Block diagram of DPG. 

The boundary conditions of the five fuel nodes required five 
values of K(T) . However only four DFG's were available on the analog 
ccmputer, consequently, one of the DFG's had to be time shared. The 
time sharing was accortplished by alternating the information to a 
DFG from two boundaries. The switching was acconplished by two 
digital/analog switches. The logic to operate one of the switches 
was supplied by a 10 kHz square wave viiich for clarity will be called 
logic variable L, ; the other switch was operated by L. the inversion 
of Lj^. The output of the DFG was fed to two track-hold arrplifiers, one 
being controlled by logic and the other by With the aid of Fig. 

4b it can be seen that viiile one switch is closed one anplifier will 
be in the track mode following the output of the DFG while the other 
is in the hold mode. When the logic signal inverts the switches 
reverse conditions and the amplifiers reverse modes. 

The analog corputer used for this model, CCMCOR-5000, is equipped 
with digital logic and in addition conpletely controllable by the 
general purpose digital ccmputer, SDS-9300. However it was necessary 
to develop a logic network which could control the compute mode of the 
analog ccnputer automatically using the digital ccitputer for varying 
periods of time. Furthermore at times it was necessary to bring the 
model up to operating temperature before the external forcing function 
was applied to the reactor kinetics. To realize this flexibility a 
complete logic network was developed which is shown in Fig. 4c. 
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B. DESCRIPTION OF DIGITAL COMPUTER SOFTWARE 



In the previous section it was pointed out that the solution to 
the heat flow for the canplete fuel pin involved one set of equations 
by applying the analog circuitry of Fig. 4a on a time sharing basis. 
The digital computer with its logic capabilities can control this time 
sharing very effectively. Essentially what occurs is that the analog 
computer is used as a large subroutine and V\^en called produces the 
solution for the axial segment under consideration. 

In this simulation the solution to one axial segment is treated 
on a continuous time basis during any transient condition. In order 
to treat each segment in this manner the time behavior of the inlet 
coolant tonperature must be incorporated within the solution for that 
segment. The digital computer with its memory capabilities can ac- 
complish this effectively. The digital computer through an analog- 
digital converter samples at a fast rate (250 Hz) the outlet tempera- 
tures of the segment being considered. This information is stored as 
an array within the digital computer and used as input data in the 
form of a analog voltage through a first order digital to analog 
converter during the solution of the next segment. Thus the solution 
can continue sequentially froti one axial segment to the next. 

The feedback shewn in equation (6) is an algebraic function and 
for this reason it is more convenient to have the digital cemputer 
handle this aspect of the simulation. However, it is necessary at 
this point to note that the feedback mechanism is a function of the 
average temperature. This situation requires that the average tem- 
perature be known at a given time in order that the digital cemputer 
can evaluate the feedback. However the average tatiperature of the 
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core is not known until all the segments have been considered. To 
obtain the average temperature, the average fuel and coolant ten- 
peratures of each segment are weighted to approximate the average 
temperature. Each segment is then evaluated on a continuous time 
basis with the digital computer supplying the closed-loop feedback. 
The operation is illustrated in Fig. 7 where the closed-loop syston 
for each axial segment is shown. Dioring the analysis of each segment 
the digital computer rapidly samples the logic condition, defined L3, 
on a 250 Hz signal which commences simultaneously with the signal 
which places the integrators in the compute mode (Fig. 8) . When L3 
is -6 volts the digital computer samples the change in ambient tem- 
perature for the fuel (AT^) and coolant (AT^) for the axial segment 
j being considered, stores the information, computes a feedback 
after weighting the tonperatures , and outputs the feedback just 
computed as a voltage to the analog computer (Fig. 8c) . L3 during 
the above operation has returned to 0 volts. The digital computer 
remains in a loop until L3 again has the value -6 volts whereupon 
the process in Fig. 8c is repeated. The condition of 0 volts on 
logic signal L2 occurring at time t ends the computation for that 
segment. For clarity the analysis of all ten segments in this manner 
will be referred to as iteration 1. 

After each segment is analyzed the fuel and coolant temperatures, 
each represented by a 250 x t matrix, are stored on a disc due to 
storage limitations within the digital computer. After all ten seg- 
ments have been analyzed the 10 x 250 x t array for each temperature 
is processed to cibtain a 250 x t array representing the average tem- 
perature as a function of time. A K^(t) is then computed based on 
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That portion enclosed by dashed lines indicates the part 
being haindled by the digital computer. 



Figure 7. Closed-Loop Feedback Schoretic for 
Each Axial Segment. 
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this average tanperature. A new iteration is made using the K^(t) 
just computed. The changes in fuel and coolant temperatures are 
sampled as before. However this iteration differs from iteration 1 
in that the loop is not closed but open. Figure 9 shows the sche- 
matic for the open- loop process. At the conclusion of analyzing the 
ten axial segments a new (t) is computed based on the temperature 
profiles just developed. The open loop process is then repeated 
until convergence of the neutron density N(t) occurs. 

When applying the open-loop iterative method, a study should be 
conducted to insure that convergence will occur before attempting 
this procedure. Proof that this system will converge is given in 
Ref. 5 and 18, however Section V of this report covers a discussion 
on the rate of convergence. 

Figure 10 represents a flow chart for the digital software program 
just outlined. 

The method for obtaining the weighting factors for the coolant and 
fuel tanperatures is now described. Normally to start the probl^ 
from a known steady-state spacial power distribution, the initial con- 
ditions on the itegrators were needed. These were easily obtained by 
putting in a step input of heat and allowing each segment to cone to 
steady state (operating level) . The final voltages on all integrators 
were sampled and the info 2 cnation put on paper tape for future use as 
initial conditions. After each segment had reached operating level, 
it was disturbed with K^(t) as a ramp input for a short period of time 

.1 sec) without feedback. The changes in average fuel and coolant 

th 

tanperatures for this time interval were sampled for each j axial 
segment and frctn these tertperatures the changes in average coolant 
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(AT^) and fuel (AT^) temperatures for the fuel pin v?ere calculated. 
The weighting factors, and W^, were calculated using the following 
relationships , 





ATc 



These weighting factors were then treated as constants for the 
duration of the transient condition. 
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IV. DYNAMIC CHECK OF MODEL 



When a physical system is described by a model it is important 
to know whether this model is adequately describing the system. 

To ascertain this, parts of the mathematical model were ccnpared 
with known mathanatical solutions or with digital ccanputer solutions. 

The approximated neutron kinetics equations (5) were ccitpared 
with the exact solution for the case of the one delayed onitter model 
using a step input of reactivity of magnitude 50 cents. Representing 
this step input by the solution of the reactor kinetics can readily 
be obtained for this input condition [Ref. 19] and is given in equation 
(22) for the one delayed emitter model. 



When practical numbers are used in equation (22) product terms con- 
taining X can be neglected and the solution can be expressed as 



sec., therefore it is unnoticable on any time scale of interest and to 




( 22 ) 



where 




c- 




(22a) 



Since 




the second exponent has a time constant 2 x 10 



-4 



44 



a good approximation it can be neglected. What occurs with this 



delaying .001 sec. Figure 11 shows both the analog solution of 
equation (5) and the exact solution given in equation (22a) . From 
the results it can b© seen that the approximation used to develop 
equation (5) appears valid. 

Equation (5) was also solved on the analog conputer using a 
triangular pulse input for K(t) as follows: 



A digital conputer integration solution utilizing the Runga-Kutta 
method with a Adams-Moulton predictor corrector with error check 



also used to solve the reactor kinetics without the approximation, 
equation (3) . The results of both solutions are shewn in Fig. 12. 
Table I shows a cenparison between the values for the two solutions 
for every tenth of a second. The agreament is quite good. 

The heat flew equations were solved on the analog conputer for 
a homogeneous material with constant thermal properties. The exterior 
of the material was kept at a constant temperature and a constant heat 
source was applied within the material. The solution for this situ- 
ation, recalling that the equations apply for a circular cylinder, is 
[Ref. 22] 



assumption is that N (t) is a step jump to | _ ^<^ instead of 




. S’i cloUovs^ 0 1 t i 1 sec- 



-4 

(absolute error < 10 per integration step) [Ref. 20 and 21] weis 
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Figure 11. Cemparison of Analog And Exact Solution For Neutron 
Kinetics For a Step Input of 50 Cents. 
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Figure 12. Catparison of Analog and Digital Solution of Neutron 
Kinetics for a Triangular Pulse Input for K (t) . 
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TABLE I 



COMPARISON OF DIGITAL AND ANALOG SOLUTION OF REACTOR KINETICS 



Time 


Digital solution 


Analog solution 


0.1 


1.054 


1.052 


0.2 


1.118 


1.117 


0.3 


1.194 


1.192 


0.4 


1.284 


1.282 


0.5 


1.393 


1.390 


0.6 


1.525 


1.522 


0.7 


1.689 


1.685 


0.8 


1.895 


1.891 


0.9 


2.157 


2.150 


1.0 


2.501 


2.481 


1.1 


2.398 


2.376 


1.2 


2.295 


2.275 


1.3 


2.194 


2.174 


1.4 


2.095 


2.075 


1.5 


1.999 


1.980 


1.6 


1.906 


1.890 


1.7 


1.816 


1.800 


1.8 


1.730 


1.712 


1.9 


1.646 


1.629 


2.0 


1.566 


1.558 


2.1 


1.566 


1.558 


2.2 


1.566 


1.558 


2.3 


1.566 


1.558 


2.4 


1.566 


1.558 


2.5 


1.566 


1.558 
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where, 



q = heat applied in Watts/cm 

a = radius of cylinder 

K = thermal conductivity of material 

a = n^ root of Bessel function of order zero, 
n 

Figure 13 shows a comparison between both solutions. It can be seen 
that the solution is quite good for nodes 1 through 4 but node 5 is 
in error by as much as 25% in early values of time and 10% in final 
time. The situation is generally recognized [Ref. 23] as being a 
problan near abrupt boundary conditions. Since abrupt changes are 
not envisioned as part of any study with this model it can reasonably 
be expected this error to be less. Furthermore the average tonper- 
ature is the important quantity of interest; thus the error will be 
further reduced in magnitude. 

Prior to making any ccanputations with this model it is wise to 
carry out a static and dynamic check of the model to insure all 
components are functioning properly. 
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Figure 13. Catparison of Analog and Exact Solution for the Five Tenperature 
Nodes. 
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V. RESULTS 



The principle objective of this report was to simulate the U1EBR 
on the hybrid ccmputer and to improve on the rate of convergence in 
the open- loop iterative method. In the developnent of the technique 
used for irrproving the rate of convergence it was discovered that a 
problem related to the stability of the converged solution existed. 
Because of this problem of stability of convergence (to be covered 
later in this section) the solution of the problan is limited to a 
maximum of 2 seconds. 

Figure 14 and 15 show the N(t) response as a result of an 
external distxjrbance of one dollar per second viiich is terminated 
at the end of one second. The disturbance is intended to simulate 
the withdrawal of a control rod. The disturbance was started fran 
a steady state operating pcwer level of 

S(?) = (lOOO Siw 

L represents the length of the fuel rod vbich is 100 cm; z is the 
axial position variable. Figure 14 illustrates the process of 
convergence of the open- loop iterative technique. Iteration number 
1 was a result of a guess for the feedback. Ten iterations were 
required to obtain the converged solution. Figure 15 shews the 
response of N(t) produced after the iteration vbich uses the weighting 
factors, and this compares favorably with the converged solution 
shewn on the same graph. 

During the development of this model it became questionable 
whether the open- loop iteration schome would produce a converged 
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Figure 14. Illustration of Convergence Process. 
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Figure 15. Ccxtparisan of Converged N(t) With N(t) Developed from 
Wei^ting Factors. 
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solution. Even though there exists mathenatical proof that con- 
vergence occurs [Ref. 7 and 18] and published results [Ref. 7] 
showed the method successful, early results with this model re- 
vealed a divergent tendency for values of time greater than 2 
seconds . 

In an effort to isolate the problem a simple integrator circuit 
(Fig. 16a) was considered, since it can be shown mathematically that 
convergence is assured when using the open-loop iterative process 
[Ref. 18] . For the open-loop iterative process, the circuit in 
Fig. 16a can be schematically represented by that in Fig. 16b. The 
procedure followed was to study the behavior of the open- loop iter- 
ative process with representing the closed-loop analog feedback 

( 2 ) 

for the first iteration. The output Y was sampled by the digital 
computer and used as input for the succeeding iteration. This process 
was continued and the following observations made: 

a) Starting with a "good solution" for the behavior of E^ 

with successive iterations was as shown in Fig. 17 (greatly exaggerated) 

( 2 ) 

The region v^iere Y began to diverge from the converged solution was 
v^en the output approached steady state condition, or when the sum 
of X + Y^^--'0. The output would eventually converge but for times 
greater than 4 seconds (time constant for this example was 1 second) 
it tended more to oscillate with approximately equal error above and 
below the converged solution. 

b) For times less than four seconds a converged solution occurred 
but it could not be maintained indefinitely. Periodically (with no 
apparent consistency) the converged solution would be lost. Figure 18 
shows a converged solution and a succeeding one to illustrate this point 
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Figure 16. Diagram of Integrator with Feedback. 
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fiHf>LOCy SOLUTION/ 







Time Ciec.) 



Figure 17. Illustration of Successive Iterations for Fig. 16 
Using Open-Loop Method. 




Figure 18. Illustration of Stability Problem Associated with 
Open-Loop Method. 
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The system would undergo a convergence process to again obtain a 
solution. Periodically (again with no apparent consistency) v^le 
in the process of converging it would lose the tendency to converge 
from one iteration to the next and the process would again have to 
repeat itself. For times up to about 3 seconds, \dien X + was 
large, the output appeared unaffected and maintained a converged 
solution for repeated iterations. 

c) The rate of convergence was slow for interval of times greater 
than 4 seconds. 

To study the behavior of the neutron kinetics equations using 
the open- loop method a simple lumped-parameter model of a nuclear 
reactor was developed. The mathematical model for the reactor is 
given by equation (25) . 





(a) 


- bU)] 


(b) 


t = aT+ AN(f) 


(C) 


KU)= 


(d) 







where, T = average temperature, 

AT = change in average temperature above T^, the initial temperature, 

1/a 5 time constant of fuel. When heat is removed the temperature 
of fuel will fall to 1/e of its steady state value in time 
equal to the time constant, 

A = factor converting neutron density to pcwer. 

The Doppler effect was the only feedback mechanism considered and it was 

approximated by equation (25e) ; a closed-loop analog solution could then 

be obtained. Again the above equations can be schematically represented 
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as Fig. 16b. The open- loop iteration solution was attenpted using 
the analog solution for Y . Similar behavior to the observations 
made above was noticed. 

Two factors to be considered when using the open-loop method 

based on the above observations are the rate of convergence and the 

stability of the converged solution. The rate of convergence for the 

systan described in Fig. 16 is slew. Ihis can be shewn as follews. 

If X (Fig. 16b) is considered as a step inpat of magnitude one, then 

for an initial guess of (t) = 0 the output, Y^^^ (t) , is equal to 

Gt (recalling the gain of the integrator is G) . For a unit feedback 

2 2 

the input for iteration 2 is 1 - Gt, and the output is Gt j-. It 

is easy to show that a continued iteration process will develop the 

follewing series, 
at/) 

= Z 3! 

vhere i denotes the iteration. A necessary condition for the series 

r 

to converge is lim 9 . q . As long as G and t are finite this con- 

i-Hx) I '■ 

dition is nvet. However before convergence begins to occur the ratio 

th 

of the i + 1 term to the i^' term must be less than one, or 

<1 (26) 

The value of G used in Fig. 16 was 10. Thus convergence can be seen 
as a slew process. It is also possible to understand the initial 
divergent tendency of the open-loop method oatmented on earlier. 
Recalling that Y^^^ in Fig. 17 represents the closed- loop systan 
solution used for the first iteration, then if their exists any dif- 
ference between the step input X of Fig. 16b and the feedback Y^^^ 
near steady state the integrator will continue to show a change in the 
output with the later values of time associated with the larger 
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differences. If it is assumed that any difference bet^/een the analog 
solution and the digital solution is constcint in time, it can be seen 
fran the criteria used in developing equation (26) that the solution 
for early times will begin to converge while later values will tend 
to diverge. Eventually all times will begin to show convergence. 

This is illustrated by iterations 1 and 3 in Fig. 14. 

There obviously exists a stability for foe converged solution in 
the open-loop method. Since the behavior of this stability problen 
was not consistant nothing specific could be found that attributed 
to the problem. The digital to analog output was adjusted to ensure 
no bias voltage existed. The error between the digital and analog 
solution for was conpared and the difference was within the re- 
solution capabilities of the equipment. The differencing of X ard 

was acconplished within the digital cotputer to nullify any error 
in the analog to digital circuit. The sampling rate at the analog to 
digital interface was varied from 250 to 400 Hz (400 Hz was the maximum 
attainable) . In each case the behavior reported for Fig. 16 was 
observed . 

The stability problan appears more pronounced when 
The neutron kinetics, especially for a fast reactor with 
are very sensitive to K(t) and if this quantity does not go to zero in 
the limit of steady state then the equations respond accordingly. It 
is a necessary pxrocedure in analog simulation to ensure that with 
K(t) =0 the reactor kinetics remain in the steady state condition vben 
in the compute mode. With equation 3 which requires an integrator with 
a gain of 10^ the system requires constant checking to insure a 

balanced condition between N and D so that the steady state solution 
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does rjot drift in the oarpute mode. The use of equation (5) , vAiich 
ranoves the integrator with gain eliminates this need for constant 
checking, and it was hoped that the stability and convergence problem 
would be reduced. Howsver, no noticable im p ro v ement occurred. 

It can be seen from Fig. 15 that the algorithm developed in this 
model rapidly produces N(t) especially for times less than 1 second. 
However, an error exists for times greater than 1 second between the 
two curves and more iterations should be carried out to reach the 
converged curve. A converged solution is obtained in agreement with 
that obtained in Fig. 14, however, an unusually large number of 
iterations are required \dien catpared with the rate observed in Fig. 

14 in the same time region. Reference 26 reports and demonstrates 
that the rate of convergence is enhanced v^en starting from a "good 
initial" guess than from a "poor" one. It is suspected that the 
stability problem discussed earlier is effecting the rate of conver- 
gence for times greater than 1 sec. More research on this prdblon 
is needed. 
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VI. CONCLUSIONS AND REMARKS 



The open- loop iteration process as a solution to partial dif- 
ferential equations is a relatively new algorithm for application 
on a hybrid computer. Presently only two reports [Ref. 7 and 24] 
have been located which present investigations into the algorithm. 

It is interesting to note that Ref. 24, which used the open-loop 
method to solve the one dimensional wave equation and heat flow 
equation, addressed itself to results short of the steady state 
solution. Also in Ref. 8, which was a follow up study on Ref. 7, 
the open-loop method was discarded because a large number of iter- 
ations were required near steady state conditions, and a closed- loop 
iteration process was adopted. 

For the open-loop iteration process to be a successful tool in 
solving partial differential equations two problems need further in- 
vestigation, namely, the rate of convergence and the stability problem 
observed in this report. The neutron equations can be reduced to an 
integral operator [Ref. 7] and therefore the convergence is slow 
especially for times near the steady state solution. A large number 
of iterations detract fron the advantages of using a hybrid computer. 
In Ref. 24 the number of iterations in viiich convergence occurred 
over a given time interval increased si±>stantially as the solution 
approached steady state. 

The problem of stabilii^ near steady state solution appears to be 
attributed to the algorithm or the hardware. It is doubtful if the 
nodel is producing the problem since similar behavior was observed in 
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a very simple system. Fig. 16. Further investigation is needed to 
establish the source of this problan. It is suggested that based 
on the investigation carried out in this report that the problem 
is inherent with the open-loop method. 

One possibility looks encouraging for future studies in solving 
the rate of convergence and the stability problem. The technique 
QTiployed in Ref. 8 utilizes a closed-loop iteration process. The 
method rapidly develops N(t) for later values of time but N(t) for 
intermediate values of time required more iterations. The algorithm 
in this report develops N(t) for early times but requires additional 
iterations to reach a solution. A ccmbination of the two may prove 
fruitful . 

In addition the use of the Diode Function Generator represents 
a very flexible tool when used to duplicate a nonlinear function. 
This is especially demonstrated when incorporated in the continuous 
solution to the nonlinear heat flow equation. The procedure used in 
Ref. 7 was to approximate K(T) as a third order polyncmial which 
resulted in the need for a large number of multipliers and summer 
amplifiers. Also the DFG can approximate the curve more accurately 
then can the third order polynomial. 

It is also worthwhile to note that the present hybrid simulation 
has the following limitations: 

1) The bonding between the fuel and cladding was neglected. 

This will produce lower temperatures in the fuel than actually en- 
countered. An attempt was made to include the bonding but the 
analog circuit was unstable. 

2) The point model neutron kinetics are used to describe the 
behavior of the entire core. 
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3) Axial flew of heat is neglected. This will result in 
producing higher temperatures. However any error in temperature 
introduced here is considered negligible. 

4) The assumption that an average fuel cell gives all the 
information necessary to obtain the average core tonperature may 
lead to unrealistic dynamic behavior [Ref. 6] . However in the 
design of the General Electric LMFBR the fissile fuel is varied 
radially within the core to flatten the power flux. This would 
tend to make the assumption more realistic. 

Presently this model is very complex and some problems with the 
present model need further investigation before the limitations 
outlined above can be considered. 
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APPENDIX A 



A. MATHEMATICAL PROOFS AND DERIVATIONS 
1. Kinetic Equations 

The kinetic equations are recast here in the reduced variable 
form. Equation (1) is considered with one precursor concentration, 
namely, 

rt . /.\ , r . . . 
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Equation (la) is divided by lf\{0) and equation (lb) by c (0) , 

3 

vdiich are the values of neutron and precursor concentration (#/c3n ) 

n ( t) 

respectively at steady state. Defining two quantities, N(t) = 



and D(t) = 
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equation (1) becomes 
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(27) 
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At steady state ~ — k ~ Q in equation (1) , therefore, 

^/i V}(h) = \ C.(a'). Putting this into equation (27) the follcwing 
equations are obtained. 
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The quantity = K(t) has the units of dollars. Thus one dollars 

worth of reactivity indicates p (the excess reactivity) has the value 
of 6 . Introducing the term K(t) equation (2) is obtained. 



L [}<»]N(i)-Nm f DU)'j (a) 

X [NM-bCi}] (W 

v^ere at t = 0, N(0) = D(0) =1. 

2. Reactivity Feedback 
a. Doppler Effect 

In a IMEBR the Doppler effect contributes the majority of 
the feedback [Ref. 9] . A highly simplified picture [Ref. 25] of the 
Doppler effect is illustrated in Fig. 19. The cross-section curve 
v^ich shows the probability of a nucleus-neutron reaction is shewn 
at left. At the right (Case I) is a picture of a neutron approaching 
a nucleus at a low tanperature. Since the nucleus is essentially at 
rest neutrons with a lesser or greater velocity than that occurring 
at resonance will escape the reaction. However as the temperature 
increases (Case II) the nucleus is no longer at rest and neutrons of 
Icwer or higher velocity can still achieve the same relative velocity. 
The effect is to broaden the resonance peak and increase the chances 
of a neutron being absorbed. The absorption of neutrons for other than 
fission processes deminates and a negative reactivity occurs. 



65 




A/eutrovi 

Veioc-ity • ^ 

Relative ^ 

Vcloci+y • 

CASB T 

• 






hJaclcui 

o 



<-o 

c>-^ 



CASE U 



Figure 19. Illustration of Doj^ser Effect. 
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A development of an exact imthematical model for the 
Doppler effect has been hampered due to lack of experimental data 
for neutron cross sections at higher energies. It is necessary to 
use statistical treatments [Ref. 26] v^ich yield for the Doppler 
effect the expression: 

dT 

However recent results by R. Froelich, and others [Ref. 
27] indicate the Doppler effect to be 



cTT 

vdiere C and y are positive constants depending on volume ratio of 
UO^ to PO 2 and sodium density, y is approximately equal to one. 

In the General Electric reactor under design y is equal to one. 
Therefore the Doppler effect is: 

t'T "t” 

vdiere T is the average temperature of the core and is a position 

constant termed the Doppler coefficient. Starting at an average tem- 
perature of the core, T , the change in reactivity due to Doppler 

o 

effect (Pj^pp) t)e e>q>ressed as follows: 

Pbopp = “ J Y ' ~ ( 1+ 

or, in terms of dollars 

b. Fuel Expansion 

The expansion of the fuel gives rise to a negative feed- 
back. It can be explained scmes^^at simply by observing that an 
expansion of the fuel reduces the fission nuclei per cm seen by the 
flux of neutrons. In the General Electric design LMFBR it is expressed 
as 
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or, in dollar units 








vdiere B is a positive constant, 
exp 

c. Coolant Expansion 

The expansion of the coolant can give rise to either a 
positive or negative feedback. As the coolant expands the nuniber of 
neutrons captured parasitically is reduced thereby increasing reac- 
tivity. Changes occur in neutron leakage (loss of scattering events 
by sodium) v^en the coolant expands and, hence, reduction of reac- 
tivity occurs. Also less moderation occurs producing a hardening 
of the neutron spectrum, and if fertile materials are present more 
fissions occur producing an increase in reactivity. The si 2 e of the 
core, fuel and volume ratio of coolant determine which of the above 
are more pronounced. In the General Electric reactor the feedback 
due to coolant expansion is positive and expressed as follows, 



\diere is a positive constant usually referred to as the sodium 

void coefficient. 

3. Heat Flew Equations 

Fran Section IIA the equation for heat flew in the fuel rod 




or, in dollar units. 



K, 





( 9 ) 
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Figure 20 illustrates the various radial segments selected within the 
fuel pin. The tenperature between mid-points of each radial segment 
is assumed to vary linearly. The slope at each side of a segment 
interface are equal thus ensuring that the flow of heat across the 
boundary is continuous. Applying equation (8) to the region indi- 
cated in Fig. 20 and integrating over the special variable r the 
following equation is obtained,^ 



'oil / 

1 n-f 

For clarity each integral will be handled separately. Assuming K(T) 

is not a constant, the first expression on the right side becatvss, 

1 1 * 















(29) 






where, = thermal conductivity at the boundary between 

node i + 1 and i. The tanperature at the boundary is -■ 

A 

The second integral on the right side reduces to. 



Reference 15 states this technique of reducing a partial dif- 
ferential equation to an ordinary differention equation requires a 
minimum number of nodes. 
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Figure 20. Teannperature Profile Within Fuel. 



Tr 




Figxare 21. Tenperature Profile Within Cladding and Coolant. 
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(30) 



V ' 






V, i i>l 

' X. 






2 



V, - 



r - C' 



^ SC?; Aya) V; .'/ 



The quantity in the integral on the left hand side of 

g 

equation ( 9 ) is assumed to be a constant , therefore/ c^d^ and the 
operator t can be brought outside the integral. The 

integral is then handled as follows: 



r, ■# 

/ aV - 
/ ^ .)t 




•a 



(30 



Fran Fig. 20 the relation for T is. 




( /,-e^ )f /c 
V. < (Yc + ""i) 



Therefore equation (31) becones. 



Ci- c// 



J 

It 




^ ( 



'f-r: 







dr 




In most materials this can be considered a constant and in UO» 
the variation is not pronounced. Reference 15 treats the condition 
vdiere c^d^ is not a constant. 
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After applying limits and seme reduction equation (32) becemes: 



dt 

If I {s 



2^ J 



(33) 



^ r 
\r- — 

3 " 2f 



\ / 

The operator in equation (33) can be changed to since the 

space variable has been eliminated. Combining equations (29) , (30) , 

and (33) and after dividing through by the following equation 

is obtained. 



4 i [t.. (i- 

- Kf... ( - o: 




(34) 



■h SU) N(i) 



The equation for heat flew in the cladding is, 
/ ^ 






( 10 ) 
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For the cladding which for this model is stainless steel the 



quantities c d and K are assumed constant. Foliating an identical 
CO o 

procedure as above and utilizing Fig. 21 equation (10) can be reduced 
to, 



To differentiate one axial segment fron the next, superscript 
j is used and equations (34) and (35) becone equations (11) and (12) 
used in Section IIB. 

To the clad-coolant interface the conservation of energy is 
applied to develop a differential equation which couples each axial 
segment. As stated earlier the axial flew of heat is neglected and 
therefore no coupling of one axial segment with its neighbor exists 
within the fuel and cladding regions. The coupling occurs through the 
flcwing coolant. 

Fig^Jre 22 represents a top and side view of an axial segment. 
If the coolant is treated as a heat sink, the temperature drops at 
the clad-coolant interface to a value assumed constant (Fig. 21) over 
the region which receives heat from the segment in question. 
represents the extent in the radially direction of this region. The 
radial gradient of tonperature in coolant is neglected. 




(35) 



The heat (6q) received by the coolant is 



bat, dTUt.t)= 'Hfdt 

(){ 



(36) 
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Top View 




Side View 



Figure 22. Diagram of Axial Segment. 
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L'OO^^A/ T 



therefore equation (36) beccttves 



cf 



1 






= 4c. [ 



c)i 




(37) 



or. 



±L 



at d~t 



JT. 

Ji 



d \ \ 

Ji J_^ 

<^Coa, / 



TT (38) 



^ 1 ^ = Mass flow rate of coolant = G(t) 

cool dt 



h 



ZTTr(_a ^ dt 



= Heat flux across clad-coolant interface h h 



Both numerator and demtdnator of the left side of equation (39) 
are nov multiplied by 2Ttr^ in order to cast it in a fom to take ad- 
vantage of certain design criteria. In a reactor the volxme fraction 
of coolant is a known design criteria and this knowledge allows the 
evaluation of 27rr^ ^^cxx>l' approximately the cross section 

of coolant. In most IMFBR designs [Ref. 28] considered to date, includ- 
ing the General Electric model used herein, the volume fraction of 
coolant ^i:50% iirplying the cross section of coolant equals the coribined 
cross section of fuel and clad. Therefore, 



ITT 

- = ratio of circumference of clad to cross sectional 

* area of coolant 



TT ~ 

Equation (38) may now be written as. 



^Jc ^ ()Tc Z H 

^Coo, dt 



(39) 
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Again an ordinary differential equation is sought in order to 



use the analog oanputer. To reduce equation (39) to an ordinary dif- 



ferential equation the finite differencing technique [Ref. 29] is 

erployed by discretizing the space variable. The operator is 

-T"'' 

reduced to — — — . Also, a linear relationship is assumed for 

T across the Az segment, thus, the average coolant taiperature (T ) 
c c 

for a segment can be expressed as 



j = 1, 2, 3 10 (10 being the number of segments) 

Ten axial segments were found adequate to describe the coolant ten- 



perature profile in agreement with Ref. 7. 



nj-1 



T^ “ represents the outlet 

torperature fron the previous segment and T^ the outlet temperature of 

o 

the segment in question. 

With these substitutions equation (40) becones. 



ciji cl^oof 



iV-V'j ZH 



^ 2 



^ o 0 > ^Coof ^ 



but, t 3 = IT‘-X 



J-l 



therefore , 



Z GU) 

dt ^ 2 cl, 



Coo ( 



( Tj - T;'-' j = 



‘^Caol ^<»el 



(14) 
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APPENDIX B 



A. SCALING OF DIFFERENTIAL EQUATIONS 

The scaling of the differential equations (15) through (21) 
will be daivonstrated by using equation (15) 

+.24.^3"h' J= (T2!'V)]-t.^looSC?)N(t) (15) 

The various quantities (T^ , T^, represented on the 

analog oorputer as voltages v^ereas in equation (15) they represent 
units (°K, watts/an^, etc.), therefore a coefficient called a scaling 
factor is needed which transforms the various units into voltages. 
Thus the following expression for each variable in equation (1) is 
defined. 

T = 

N(^}= k, nU] 

s^?)= ks sCi) 

vihere T^, T 2 , ^ volts. Placing these 

quantities into equation (15) and dividing through by gives 

i [T,’+ • ]. 3is.z [k. (Tj- T')] 

(40) 

f .^1 ) s(?) 

^ lej / 

The value of the scaling factors k^, k^^, and k^ are determined by 
knowing the maximum voltages the variables will assume during the 



77 



solution. Ihese are sometimes not known and often the problem must 
be put on the ocnputer and tested before arriving at scaling factors. 
In this case it was found convenient to use the following values for 
the scaling factors, 



IJt= 2d'' ^ Aon 






= 



Watts 



OO 



Volt 



- * Vvoit 

With these scaling factors equation (40) beccnies 



+. 2t3J li 3Z,8z[/o ~ sCi)N(i) ( 41 ) 

Similar operation on the remaining equations (16) through (21) produce 
the following equations 

[1;'+.Z=38T,' I21b-T']. 

. _ . - 

- (V-7‘)j -I- .2 Ui- fa) NU) 



[ 1-4 + • + . /444T,' J = lo.S [/o (tJ-TJ)] 

- [lo k^3 2 (~i~A~^i)]+ . 2iti" s^:?) nU) 

JiiV /<rz4Tj ^,/So 7T/]= ic,,ai !--%')] 

- /4.8b' [/o (%‘-V)] + .llbS- l(z) NCt) 



(43) 





1^0 Ih 



— /^>7i‘[/o /< 






(4S) 
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- &<?.s £/0 /?^,fc (%'-Ts^)] 

T^ --^ 5 - 4 . 3 Tt •^ IIS.sTg -I- 337. s ' (-^7^ 

(Note: T ^ refers to average coolant tonperature in volts) 

W 
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APPENDIX C 



(7) 

A. PHYSICAL PARAMBIERS ^ ^ 

Specific heat-c (watt-sec/gm-°C 
Fuel - 0.3345 
Cladding- 0.502 
Coolant - 0.21 

Thermal oonductivity-K (T) (watts/an-°C 
Fuel - see Fig. 5 
Cladding - 0.2075 
Density-p (grry/an^) 

Fuel - 9.2 
Cladding - 8.0 
Coolant -0.8 

2 

Mass flav rate of coolant-G(t) (giVon -sec 
500 

Heat transfer coefficient between clad and coolant-h (watts/an^-'^C 

13.08 

Dimensions (an) 

Ar = 0.05 

Ar = 0.05 
c 

Az = 10 
L = 100 

ProTpt neutron lifetime-^ (sec) 

4.8 X lO"^ 



Values taken frcm Ref. 9 and personal correspondence with Mr. C. K. 
Sanathanan. 
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Decay cx>nstant for six precursor nodel-Xi(sec 

L — ^ 

1 0.0130 

2 0.0314 

3 0.136 

4 0.340 

5 1.320 

6 3.500 
Delayed neutron fraction from group-gi 



1 

1 

2 

3 

4 

5 

6 



.0000759 

.000626 

.000564 

.0010700 

.000489 

.000163 



& = ^ f>i = .002988 

' L-l 



- 1 , 



Decay constant for one precursor model- X (sec ) 



A 



71 ^ 

TJT 



= .5808 






= 0.005 



B 



exp 

'cool 



= 3 X 10 



-7 



= 2 X 10 



-6 
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